#include<iostream>
using namespace std;
#include<cmath>
#include<fstream>

#include"library.h"

void derivative(int n, double x, double *y, double *dydx){
	dydx[0] = 1;// dx/dt = v = 1
	dydx[1] = y[2]; // dy/dt = y^{(1)}
	dydx[2] = -10; // dy^{(1)}/dt = -10
}

int main(){
	int i,j, n=3, N=10000;
	double a=0, b=10;
	double h = (b-a)/N;
	cout<<"h="<<h<<endl;
	double *y = new double [n];
	double *ynext = new double [n];
	ofstream fp("exe2.txt");
	y[0]=0; y[1] = 0; y[2] = 0;
	fp<<a<<"\t"<<y[0]<<"\t"<<y[1]<<"\t"<<y[2]<<"\n";
	for(i=0;i<N;i++){
		rk4(a+i*h, h, n, y, ynext, derivative);
		for(j=0;j<n;j++) y[j] = ynext[j];
		fp<<a+(i+1)*h<<"\t"<<ynext[0]<<"\t"<<ynext[1]<<"\t"<<ynext[2]<<endl;
		if( fabs( y[1] + 5 ) < 1E-6 ){
			cout<<"When y = -5m, x = "<<y[0]<<"m.\n";
		}
	}
	fp.close();
	cout<<"When t = 10s, x = "<<y[0]<<"m, y="<<y[1]<<"m.\n";

	delete [] y; delete [] ynext;

	return 0;
}
